Decision Trees

O Problema

Esse notebook busca introduzir o algoritmo de Árvore de Decisão, para isso será usado o dataset Iris.

# estabelecendo o ambiente
suppressMessages({
        library(tidyverse)
        library(magrittr)
        library(knitr)
        library(GGally)
        library(rpart)
        library(rattle)
        library(rpart.utils)
        })  

setwd("~/Dropbox/kaggle/iris-species/")  
opts_chunk$set(cache=TRUE)  

data(iris)
iris %<>%  as_tibble()

O dataset contém medidas de comprimento e largura da pétala e da sépala de três espécies da flor Íris:

Versicolor

Setosa

Virgínica

Entendendo a flor:

Uma amostra do dataset, ao todo ele contém 150 amostras, 50 de cada espécie.

iris %>%  head() %>% print()

## # A tibble: 6 × 5
##   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
##          <dbl>       <dbl>        <dbl>       <dbl>  <fctr>
## 1          5.1         3.5          1.4         0.2  setosa
## 2          4.9         3.0          1.4         0.2  setosa
## 3          4.7         3.2          1.3         0.2  setosa
## 4          4.6         3.1          1.5         0.2  setosa
## 5          5.0         3.6          1.4         0.2  setosa
## 6          5.4         3.9          1.7         0.4  setosa

A missão

Conhecer o algoritmo de Árvore de Decisão e cirar um modelo que possa prever a espécie de Íris a flor a partir das medidas da pétala e sépala.

O primeiro passo é dividir nosso dataset em sets de treinamento e de teste. Para que as análises daqui pra frente sejam feitas somente no conjunto de treinamento e não influenciem as medidas feitas sobre o testset.

set.seed(654)
train_idx <- sample(nrow(iris), .75*nrow(iris))
train <- iris %>% slice(train_idx)
test <- iris %>% slice(-train_idx)

Ok, ao todo temos 4 preditores, vamos checar como os mesmos estão relacionados entre si.

train %>% ggpairs(aes(colour=Species), columns=1:4, lower="blank", diag="blank", 
                  upper=list(continuous="points")) + theme_bw()

Aprendizado de Máquina

Dentro do espaço Preditor \(\Large{X}\) encotrar a melhor subdivisão que explique as classes, nesse caso as espécies de Íris.

\[ \Large{\mathcal{F}} : X \rightarrow Y \]

Árvore de Decisão

Método de classificação e regressão não-paramétrico e não-linear.
Envolve estratificação e segmentação do espaço de preditores em um número de pequenas regiões.

Definições

Abordagem Top-Down -> imagina a árvore invertida, o cume do grafo é a Raiz.
Nós internos -> Ramais.
Nós terminais -> Folhas, cada qual designa uma classificação e são as micros regiões oriundas das divisões feitas no espaço dos preditores.

tree.fit <- rpart("Species ~ Petal.Length + Petal.Width", train, control=rpart.control(cp=0, minbucket=1))
fancyRpartPlot(tree.fit, sub="")

Regiões de decisão

limiares <- rpart.subrules.table(tree.fit) %>% 
                as_tibble() %>%
                select(Variable, Less) %>% 
                filter(!is.na(Less)) %>% 
                mutate(Less=as.numeric(as.character(Less)))
                
limiares_pl <- limiares %>% filter(Variable=="Petal.Length") %>%  .$Less           
limiares_pw <- limiares %>% filter(Variable=="Petal.Width") %>%  .$Less           

x_axis <- train$Petal.Length %>% range()
y_axis <- train$Petal.Width %>% range()
n_points <- train %>% nrow()
            
train %>% ggplot(aes(x=Petal.Length, y=Petal.Width, colour=Species)) +
            geom_point() + theme_bw() + scale_x_continuous(breaks=seq(10,1)) +
            geom_line(aes(x=rep(limiares_pl[1], n_points), y=seq(y_axis[1], y_axis[2], length.out=n_points)), colour="black") + 
            geom_line(aes(x=seq(limiares_pl[2], x_axis[2], length.out=n_points), y=rep(limiares_pw[1], n_points)), colour="black") +
            geom_line(aes(x=rep(limiares_pl[2], n_points), y=seq(y_axis[1], y_axis[2], length.out=n_points)), colour="black") +
            geom_line(aes(x=seq(limiares_pl[2], x_axis[2], length.out=n_points), y=rep(limiares_pw[1], n_points)), colour="black") +
            geom_line(aes(x=rep(limiares_pl[3], n_points), y=seq(y_axis[1], limiares_pw[1], length.out=n_points)), colour="black") +
            geom_line(aes(x=seq(limiares_pl[3], x_axis[2], length.out=n_points), y=rep(limiares_pw[2], n_points)), colour="black")

Cada região dessa representa uma folha. O espaço preditor foi segmentado aqui, porém se pusermos um exemplo em 3D como o de abaixo, vemos a estratificação.

Um exemplo aleatório:

Vista superior:

Outro detalhe importante sobre a divisão do espaço de preditores é o que não acontece, como abaixo:

Como construir a árvore?

– Ideia Básica:

Escolha o melhor preditor \(A\) e o threshold para raiz da árvore.
Separe o trainset \(S\) e, subsets \({S_1, S_2, ..., S_k}\), onde cada subset \(S_i\) contém amostras que seguem a mesma condição dada.
Recursivamente, aplique o algoritmo para cada novo subset \(S_i\) até todos os nós terminais contenham a mesma classe.

Como escolher o melhor preditor feature no nó da árvore?

Taxa do erro? -> Não converge bem, não é sensitiva ao crescimento da árvore.
Índice Gini -> Medida de pureza do nó.

\[G({s_i}) = \sum_{k=1}^{K}{\hat{p}_{{s_i}k}(1-\hat{p}_{{s_i}k})}\]

Entropia -> Quantidade de desordem.

\[E({s_i}) = -\sum_{k=1}^K{\hat{p}_{{s_i}k}\log{{\hat{p}}_{{s_i}k}}}\]

Comportamento da função em relação a probabilidade das classes:

f_gini <- function(p){ p*(1-p) + (1-p)*(1-(1-p)) }
f_entr <- function(p){ ifelse(p%in%c(0,1), 0, 
                              - (p*log(p, base=2) + (1-p)*log((1-p), base=2)))}

ps <- seq(0, 1, length.out=100)
y_gini <- sapply(ps, f_gini)
y_entr <- sapply(ps, f_entr)

ggplot(tibble(probs=ps)) +
        geom_line(aes(x=probs, y=y_gini, colour="Gini")) +
        geom_line(aes(x=probs, y=y_entr, colour="Entropia")) +
        theme_bw()

*Exemplo divisão de uma região

1 - Medição de entropia no Nó Raiz, no caso usando somente os preditores: Petal Length e Petal Width.

O espaço Preditor é a região mostrada abaixo, no caso \(S_1\).

S1 <- train %>%  select(Petal.Length, Petal.Width, Species)

ggplot(S1, aes(x=Petal.Length, Petal.Width)) + 
        geom_point(aes(colour=Species)) +
        theme_bw()

Para calcular a entropia, precisamos somente da probabilidade de cada classe.

\[E({s_i}) = -\sum_{k=1}^K{\hat{p}_{{s_i}k}\log{{\hat{p}}_{{s_i}k}}}\]

S1 %>% group_by(Species) %>% 
        summarise(quantidade=n()) %>%
        mutate(prob=quantidade/sum(quantidade)) %>% 
        mutate(prob=round(prob, 3))

## # A tibble: 3 × 3
##      Species quantidade  prob
##       <fctr>      <int> <dbl>
## 1     setosa         38 0.339
## 2 versicolor         38 0.339
## 3  virginica         36 0.321

Portanto a entropia daquele nó, pode ser calculado como:

E_S1 <- - (0.339*log(0.339, base=2) + 0.339*log(0.339, base=2) + 0.321*log(0.321, base=2)) 
print(E_S1)

## [1] 1.584349

2 - Critério para divisão:

Ganho de Informação: Soma ponderada da entropia dos nós filhos ser menor do que a do nó pai.

\[\Delta{E} = p(S_{1})E(S_{1}) - \Big[p(S_{11})E(S_{11}) + p(S_{12})E(S_{12})\Big]\]

No nosso caso temos duas variáveis contínuas. Para cada atributo, são testados todos os valores continuos do seu range. O limiar que obter o maior ganho de informação é o selecionado.

fs_entropy <- function(S){
        if( nrow(S)==0){
                return( 0)
        }
        
        list_p <- S %>% 
                group_by(Species) %>% 
                summarise(quantidade=n()) %>%                   
                mutate(prob=quantidade/sum(quantidade)) %>% 
                .$prob
        
        list_p <- list_p[list_p > 0 | list_p < 1]
        entropia <- list_p %>% 
                        sapply(X=., FUN=function(p){ -p*log(p, base=2) }) %>%
                        sum()
        
        return(entropia)
}


delta_E <- tibble(attr=character(), thrs=numeric(), gain=numeric())
E_S1 <- fs_entropy(S1)
for( A in colnames(S1)[1:2]){
        range <- S1 %>% 
                select(eval(parse(text=A))) %>% 
                range()
        range_seq <- seq(range[1], range[2], 0.01)
        
        for( i in range_seq){
                S11 <- S1 %>% filter(get(A, pos=.) >= i)
                S12 <- S1 %>% filter(!get(A, pos=.) >= i)
                
                E_S11 <- fs_entropy(S11)
                E_S12 <- fs_entropy(S12)
                
                dE <- E_S1 - (1/nrow(S1))*(nrow(S11)*E_S11 + nrow(S12)*E_S12)        
                
                delta_E <- rbind.data.frame(delta_E, tibble(A, i, dE))
        
        }
        
        
}

ggplot(delta_E, aes(x=i, y=dE)) + 
        geom_line(aes(colour=A)) +
        facet_wrap(~A, nrow=1, scales="free_x") +
        theme_bw()

Aqui vemos que existe um empate e os dois atributos atingem o valor máximo de ganho. Em petal.length esse valor de topo se repete dentro de uma faixa de valor.

delta_E %>% filter(dE==max(dE)) %>%  filter(A=="Petal.Length") %>%  select(i) %>%  range()

## [1] 1.91 3.00

O valor médio dessa faixa, é exatamente o valor de corte escolhido na nossa primeira árvore. (O número aparente no diagrama árvore é arrendado por isso não bate examente.)

(1.91+3)/2

## [1] 2.455

3 - Repetir o passo 1 e 2.

Com os Nós filhos gerados da etapa anterior, repetir o processo de divisão.

Podagem da Árvore e Overfitting

\[\sum_{m=1}^{|T|}\sum_{x_i\in{S_m}}(y_i - \hat{y}_{R_m})^2 + \alpha|T|\]

Podando a árvore modelada anteriormente.

tree.fit.pr <- prune(tree.fit, cp=0.1)
fancyRpartPlot(tree.fit.pr, sub="")

Predição

Acertividade no trainset da árvore inteira e podada.

train.pred.1 <- predict(tree.fit, train, type="class")
## acertividade no trainset
mean(train.pred.1==train$Species) %>% round(3) %>% `*`(., 100)

## [1] 99.1

train.pred.2 <- predict(tree.fit.pr, train, type="class")
## acertividade no trainset
mean(train.pred.2==train$Species) %>% round(3) %>% `*`(., 100)

## [1] 95.5

Acertividade no testset da árvore na íntegra e podada.

test.pred.1 <- predict(tree.fit, test, type="class")
## acertividade no testset
mean(test.pred.1==test$Species) %>% round(3) %>% `*`(., 100)

## [1] 94.7

test.pred.2 <- predict(tree.fit.pr, test, type="class")
## acertividade no testset
mean(test.pred.2==test$Species) %>% round(3) %>% `*`(., 100)

## [1] 94.7

Elencamentos

Construções de modelos bem mais poderosos, usando árvores.

Bagging

Bootstrap e combinação paralela.

Boosting

Combinação serial das árvores.

Random Forest

Combinação paralela de centenas de árvores com truques de decorrelação.

O Lado A e B do modelos de Árvore de Decisão

Prós

Boa interpretabilidade.
Funciona bem para todos os tipos de features (caracter factor, númerico, booleano).
Insensível a outliers.
Fácil de implementar.

Contras

Não performa bem para limites lineares ou suaves.
Tendência a overfitting (segue demasiado o ruído)
Não competitvo aos melhores algortimos de aprendizado supervisionado. Contudo existem métodos de elencamento, que melhoram significativamente sua performance ao custo de perder interpretabilidade.